Numerical Estimation of the Asymptotic Behaviour of Solid Partitions of an Integer 
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The number of solid partitions of a positive integer is an unsolved problem in combinatorial 
number theory. In this paper, solid partitions are studied numerically by the method of exact 
enumeration for integers up to 50 and by Monte Carlo simulations using Wang-Landau sampling 
method for integers up to 8000. It is shown that lim n ^oo igajswj = 1.79 ± 0.01, where pz(n) is 
the number of solid partitions of the integer n. This result strongly suggests that the MacMahon 
conjecture for solid partitions, though not exact, could still give the correct leading asymptotic 
behaviour. 
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I. INTRODUCTION 

Combinatorial enumeration problems arise naturally 
in many problems of statisticalphysics. The number of 
partitions of an integer (see Q, for an introduction) is 
one such enumeration problem with a history dating back 
to Euler. Examples of applications to physical problems 
include the q — > 00 Potts model Q , compact lattice ani- 
mals 0,3, crystal growth lattice polygons 0, Bose- 
Einstein statistics 0, Q and dimer coverings . The so- 
lution to the integer partitioning problem is known for 1- 
dimensional and 2-dimensional partitions. However, not 
much is known about higher dimensional partitions. Nu- 
merical estimation of the asymptotic behaviour of these 
higher dimensional partitions could lead to theoretical in- 
sights. In this paper, we determine numerically the lead- 
ing asymptotic behaviour of 3-dimensional partitions by 
exact enumeration and Monte Carlo techniques. 

A 1-dimensional or linear partition of an integer is a 
decomposition into a sum of positive integers in which 
the summands are ordered from largest to smallest. A 
2-dimensional or plane partition of an integer is a decom- 
position into a sum of smaller positive integers which are 
arranged on a plane. The ordering property generalises 
to the summands being non-increasing along both the 
rows and the columns. Generalisation to d-dimension is 
straightforward. Consider a d-dimensional hyper cubic 
lattice with sites labelled by i = 12, ■ ■ ■ , id), where 
if. = 1, 2, . . .. An integer height h(i) (corresponding to a 
summand) is associated with site i. A d-dimensional par- 
tition of a positive integer n is a configuration of heights 
such that 



h(i) > 0, 

h(i) > max h(ii,i2, 

Kk<d 



,U), (1) 



h(i) = n. 



The second condition in Eq. means that the heights 
h(i) are non-increasing in each of the d lattice directions. 
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FIG. 1: Partitions of the integer 4 in (a) one dimension and 
(b) two dimensions. 



As an illustration, the linear and plane partitions of 4 are 
shown in Fig. ^a) and Q^b) respectively. 

Let Pd(n) denote the number of partitions of n in d- 
dimensions. The generating function Gd(q) is then de- 
fined as 



G d (q) =J2p d (n)q n , 



(2) 



n=0 



where pd(0) = 1. The generating function for linear par- 
titions is due to Euler and is 



G 1 (q) = H(l- q k )-\ 



(3) 



k=l 



and pi (n) for large n varies as |l0| 



Pl (n)~^eKp^p j, n»l. (4) 
The corresponding formulae for plane partitions are |ll| 



G 2 (q) = 11(1- 



,k\-k 



(5) 



k=l 



and 



p 2 (n) 



C-2 



,25/36 



exp 



(a 2 n 2/3 ) , n » 1, (6) 



2 



where c 2 = 0.40099... and a 2 = 2.00945... [T^ |. 
While the generating functions for three and higher di- 
mensional partitions are not known, it is known that 
linin^oo In (pd(n)) /n d ' has a finite non-zero limit 
. We define ad to be 



ad = lim 



h\(p d {n)) 
n d/(d+i) ■ 



(7) 



In this paper, we numerically estimate a% for solid (3- 
dimcnsional) partitions to be 



a 3 = 1.79 ±0.01. 



(8) 



Generalising the results for linear and plane partitions 
to higher dimensions, MacMahon suggested that the gen- 
erating function for d-dimensional partitions could be 



TABLE I: Solid partitions for n = 26 to n = 50. 



n 


p 3 (n) 


n 


p 3 (n) 


29 


714399381 


40 


352245710866 


30 


1281403841 


41 


605538866862 


31 


2287986987 


42 


1037668522922 


32 


4067428375 


43 


1772700955975 


33 


7200210523 


44 


3019333854177 


34 


12693890803 


45 


5127694484375 


35 


22290727268 


46 


8683676638832 


36 


38993410516 


47 


14665233966068 


37 


67959010130 


48 


24700752691832 


38 


118016656268 


49 


41495176877972 


39 


204233654229 


50 


69531305679518 



Gd m \ q )=n^-^y {k - 2) - 



(9) 



fe=i 



Equation 10 is usually known as the MacMahon conjec- 
ture. Clearly, G d {q) is the correct result for d = 1,2. 

However, it is known that G d m \q) is different from Gd(q) 
for all d > 3 0,0]- In particular, in 3-dimensions 



fc\-fc(fc+l)/2 



(10) 



k=l 



gives the wrong answer for the number of solid partitions 
of 6,7,8,... . 

Let p d i(ri) denote the coefficient of q n in G B (q)- 

The asymptotic behaviour of p^ (n) for large n can be 
determined from Eq. (|10|) by the method of steepest de- 
scent. In Appendix we present a heuristic derivation 
of the large n behaviour of the coefficient of q n in the 
expansion of the infinite product 



F(q) = f[(l 



k\—a±k —a2k — a3 



(ii) 



k=l 



Substituting a\ — 1/2, a 2 = 1/2 and 03 = in Eq. I|A5J| . 
we obtain 



( m ) 1 \ 
p\ >{n) 



£3 

? 61/96 



exp 



where c' 



(m) 



(12) 



is a constant and 



4 m) 



(m) 
73 



15V43 

\/T5C(3) 



15 5 / 4 C(3) 2 



27/4^5 



1.7898..., 
= 0.3335..., 

-0.0414. 



(13) 
(14) 

(15) 



Comparing the values for a^ in Eq. © and a 3 in 
Eq. (|13fl . we conclude that the MacMahon conjecture, 



though not exact, could still give the correct leading 
asymptotic behaviour for solid partitions. The value of 
a 3 m ' 1 is a function of only a\ in Eq. i|ll|) . Thus, if we as- 
sume that the asymptotic behaviour for solid partitions 
is correctly captured by a product form as in Eq. 1)11(1 . 
then it should have the form n fe (l - (? fe )-( 1 /2±o.oi2)fc 2 _ 

The rest of the paper is organised as follows. In Scc.lTTl 
we present the results of the exact enumeration study. In 
Sec. IIIII we describe the Monte Carlo algorithm and the 
simulation results for plane and solid partitions. Finally, 
we conclude with a summary and conclusions in Sec. II VI 



II. EXACT ENUMERATION 

Previous attempts at studying solid partitions on 
the computer have been based on exact enumeration 
[l3L YwL Ho. fl7| . Tables of p^{n) exist for n up to 28 
[L5| . The table is extended up to n = 50 in this paper by 
using the standard back tracking algorithm . The al- 
gorithm is made faster by the following. Partitions that 
are related to each other by symmetry operations are 
counted only once and multiplied by the corresponding 
symmetry factor. Also, parts of the partition that are 
restricted to planes are generated by using the known 
generating functions for plane partitions. In Table [IJ we 
list the solid partitions from n = 29 to n — 50. For solid 
partitions up to n — 28, we refer to Ref. |l5| . 

We compare the exact enumeration results with the an- 
swer predicted by the MacMahon conjecture. In Fig. |3 
we show the variation of ln[p3(n + l)/p 3 (n)\ with n for 

bothp3(n) as well as p\ m (n). While there seems to be a 
good agreement, we are unable to determine the precise 
asymptotic behaviour of P3(n) from these 50 numbers. 
This is possibly due to the presence of strong corrections 
to the leading asymptotic behaviour. It is difficult to fur- 
ther extend the table of solid partitions due to the large 
computing times involved. One possible method of prob- 
ing larger values of n is to use Monte Carlo simulations. 
These are described in Sec. IIIII 
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FIG. 2: The results from exact enumeration are compared 
with P3" 1 ' obtained from the MacMahon conjecture. 



III. MONTE CARLO SIMULATION 
A. Algorithm 

We use an algorithm proposed recently by Wang and 
Landau for measuring density of states in spin systems 
|l9). The algorithm is described below. Consider a 
N x x N y x N z lattice with initial height h(i) assigned 
to each lattice point in such a way that the configura- 
tion is a valid solid partition. To each positive integer n 
is associated a histogram H (n) and the number of solid 
partitions p 3 {n). The histogram H(n) keeps track of the 
number of times solid partitions of n have been visited 
during the simulations. The algorithm is based upon the 
fact that if the probability of transition to a solid parti- 
tion n is proportional to [753 (n)] -1 , then a flat distribution 
is generated for the histogram H(n). At the start of the 
program H(n) = and P3(n) = 1. 

A site is chosen randomly and as a trial move the height 
h(i) is increased or decreased by one with equal probabil- 
ity, provided that the new state is an allowed partition. 
If the new state is an allowed partition, then the move is 
accepted with probability 



Prob(n D 



id 



>) 



Pz{n id) 
Ps{n new ) 



(16) 



where n id and n new are the sum of heights for the old 
and new states respectively, i.e, n new = n id + 1 , n i4 — 1 
or n Q id depending on whether the height increased by one, 
decreased by one or remained the same. The histogram 
H(n) and 773(71) are updated as 



H(n r - 



H(n r , 



1. 



P3 y^new ) — fiP3 (^neiu ) ? 



(17) 
(18) 



where fi is a modification factor greater than 1. 

These steps are repeated until a flat histogram is cre- 
ated; in practice, this means that H(n) m i n > cH(n) ave , 



where c is a flatness criteria typically between 0.75 
and 0.9 while H(n) min is the minimum of the H(n)'s 
and H{n) ave is the average of the ff(n)'s. When the 
histogram becomes flat, the modification factor /j is 
changed to 



t+i 



Ji ' 



(19) 



and the histogram is reset to zero. The exponent a is less 
than 1 and defines the smoothness of the iteration. The 
program runs until / is less than a predetermined value 

f final • 

Note that the algorithm does not obey detailed balance 
during the simulations. However, in the limit /, — > 1 
when p 3 (n) takes its correct value, the system does obey 
detailed balance with the weight of a solid partition of n 
being proportional to ^(n)]^ 1 . 

The algorithm can be made faster by adopting certain 
ideas from the iV-fold method [2(1, l2ll | . In this modifica- 
tion, sites which cannot undergo a valid move are never 
chosen. To do so, we define four classes, (i) cl: sites at 
which the height can only increase, (ii) c2: sites at which 
the height can either increase or decrease, (iii) c3: sites 
at which the height can only decrease, (iv) c4: an auxil- 
iary class to help implementation of detailed balance. 

First, we update the histogram H(n) and p 3 {n) by A 
times with Prob(A = k) = p' k (l —p'), where p' = ( | C4 1 + 
|d|/2 + |c 3 |/2)/(|ci| + |c 2 | + \c 3 \ + |c 4 |) and \c k \ denotes 
the number of elements in the class c k - We then choose 
one of the classes cl, c2, c3 with probabilities |ci|/(|ci| + 
2|c 2 | + |c 3 |), 2|c 2 |/(|ci| + 2|c 2 | + |c 3 |) and |c 3 |/(|ci| + 2| C2 | + 
| C3 1 ) respectively. A site is picked up randomly from the 
chosen class and a trial move is decided, for example if 
site i from class cl is chosen, the trial move is h(i) = 
h(i) + 1, while in the case of class c2 the height increases 
or decreases with equal probability. Finally, we either 
accept or reject the trial move according to Eq. i|16|) and 
update the histogram and pz(n) according to Eqs. l(T7|) 
and 118|) . With this construction, only valid trial moves 
are chosen at each time step and the algorithm becomes 
considerably faster. The role of the class c4 is to make 
the algorithm obey detailed balance asymptotically, i.e. 
when ^ -> 1. We define \c4\ = C- (|ci| + |c 2 | + |c 3 |), 
where C is some large enough constant. 

Further speeding up can be done by dividing the inter- 
val 1-8000 to smaller slightly overlapping intervals (14 
in our case) and the simulation is done for each inter- 
val. After the simulations are over, these intervals can 
be joined together to produce pz{n). The distribution is 
finally normalised by setting p 3 {V) = 1. 

The parameters we have used for the simulations are 
c = 0.85, a - 1.0/1.4, C = 2000, / = 2.5 and / /w = 
1.0000099 (corresponding to 35 iterations). Changing 
these parameters slightly does not change the final out- 
come of the simulation. The lattice sizes used for plane 
and solid partitions were 100 x 100 x 1 and 50 x 50 x 50 re- 
spectively. Random numbers were generated using stan- 
dard RANMAR algorithm. With this setup, a typical run 
producing one p 3 {n) for n between 1-8000 takes about 
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FIG. 3: The simulation results for plane partitions are com- FIG. 4: Simulation results for solid partitions are compared 
pared with the exact answer. In the inset, the variation of with p ( 3 m) obtained from the MacMahon conjecture. In the 
the relative error with n is shown. inset, we show the relative error with respect to the answer 

obtained from exact enumeration. 



12 hours with a Pentium 4 processor. For statistics we 
performed 20 runs for plane partitions and 24 runs for 
solid partitions using different random number sequences. 
Since J5 3 (n) is typically a very large number, the quantity 
that we keep track of in the simulations is ln(p 3 (n)). 



B. Simulation results for plane partitions 

We first test the algorithm against the known case of 
plane partitions. In Fig. [21 we compare the simulation 
results with the exact answer. The two curves are almost 
indistinguishable. Consider the relative error defined by 

|ln(^ s) (n))-ln(p 2 (n))| 

o{n) = 7—, -, y. , (20) 

m(p 2 (n)J 

where p^ (n) is the value obtained from simulations. We 
show the variation of 8(n) with n in the inset of Fig. [3] 
The relative error goes to zero for large n. Thus, we 
conclude that the algorithm does give the correct lead- 
ing asymptotic behaviour. In principle, the simulations 
can be made arbitrarily precise, and the correction terms 
can be determined. However, in our case, the statistical 
errors are not small enough to allow a reliable determi- 
nation of the correction terms to the leading asymptotic 
behaviour. 



TABLE II: The results obtained from the least square fit 

are shown. 0:2, c4 , a 3*™' an d a 3^ correspond to P(2,3)( n ) 
obtained from the exact results for plane partitions Eq. 0, 
Monte Carlo results for plane partitions, MacMahon conjec- 
ture for solid partitions Eq. (1101 and Monte Carlo results for 
solid partitions respectively. 





W 




(s) 


a -2 








2.010 ±0.002 


2.01 ±0.01 


1.789 ±0.002 


1.79 ±0.01 



ln[p ( 3 s \n + l)/p ( 3 s '(n)} in the range 20-8000 to the form 
0.75a 3 n- 1 /' 1 + O.b fcn- 1 / 2 + 0.257 3 tj- 3 / 4 ± d 3 n~\ where 
p 3 (n) is the value for solid partitions obtained from sim- 
ulations. We choose this form since the MacMahon con- 
jecture has the same functional form. As a test for the 
fitting routine, we test it against the simulation results as 
well as against the exact results for plane partitions us- 
ing the fitting form 0.667a 2 n- 1 / 3 ±0.333/3 2 n _2/3 ±rf2n _1 . 
The results for a( 2 3 ) are presented in Table. [H] 

For plane partitions, the value of a 2 obtained from 
the fitting routine is in excellent agreement with the ex- 
act answer 2.00945.. Agreement with the correct 

answer is also seen for a 3 m ' (see Eq. (|bHl h Hence, we 
conclude that a 3 — 1.79 ± 0.01. 



C. Simulation results for solid partitions 

For solid partitions, we calculated p 3 (n) numerically 
by averaging over 24 different runs. In Fig. 4, we 
show the results from simulation while in the inset 
of Fig. 4, we show the relative error. We estimate 
the asymptotic behaviour by fitting the data to an as- 
sumed form by the method of least squares fit. We fit 



IV. SUMMARY AND CONCLUSIONS 

In summary, we studied numerically the problem of 
solid partitions of an integer. Using exact enumeration 
methods, we extended the table of solid partitions for in- 
tegers up to 50. However, we were unable to determine 
the precise asymptotic behaviour of solid partitions from 
these 50 numbers. Solid partitions for larger values of n 
were studied using Monte Carlo simulations. From these 
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simulations, we showed that lim n —, oc> n~ 3 ^ 4 \n(p3(n)) = 
1.79 ± 0.01. This value is consistent with the MacMa- 
hon value for solid partitions. Thus, if we assume that 
the asymptotic behaviour for solid partitions is correctly 
captured by a product form as in Eq. (|llf) . then it should 
have the form T] fc (l - q k )-W±™u)k\ 
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ample, see [221), we obtain 

2oiC(4) , a 2 C(3) , a 3 C(2) 



ln(^(e- £ )) = 



e° e 
«2 + 6a 3 



12 



ln(e) + 0(e°), (A2) 



where Ci n ) is t ne Riemann zeta function. Let the coeffi- 
cient of q n in F(q) be denoted by c(n). Then 
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APPENDIX A: ASYMPTOTICS FOR THE 
MACMAHON CONJECTURE 

In this appendix, we present a heuristic derivation of 
the asymptotic behaviour of the coefficient of q n in the 
expansion of the product 



F(q) = f[(l- q 



k\—a±k —0,2 £—0,3 



(Al) 



k=l 



Let q = e~ e . Taking logarithms on both sides of Eq. (|A1|) 
and converting the resulting summation into an integral 
by using the Euler-Maclaurin summation formula (for ex- 



C(n)= 2^/^TT- 



(A3) 



For large n, we evaluate c(n) by the method of steep- 
est descent. The saddle point is the maximum of en + 
ln(F(e)). This occurs at en, where 



a i t -1/4 , \/l5a 2 C(3) _ 1/2 
en = — m-n ' H = n ' 



-3/4 



15 1 / 4 (aia 3 ^ 6 -45alC(3) ; 
24afV 



0(n _1 XA4) 



Evaluating the integral about this saddle point, we obtain 

ln[c(n)] - ^&^ 3/4 
151/43 



15a 2 C(3) 1/2 5 1 / 4 (a 1 a 3 ^ 6 -45 a iC(3) 2 ) 1/4 



/ai7r^ 



2a? /4 3 3 / 4 7r5 



11 + 21 + ^] W n ) +0 (n°). 
8 48 ° ' 17 v ' 



(A5) 
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